이 과정에서는 지리 공간 데이터 또는 지리적 위치가 있는 데이터를 다루고 시각화하는 다양한 방법에 대해 알아봅니다.
그 과정에서 다음과 같은 몇 가지 실제 문제에 대한 솔루션을 제공하게 됩니다.
글로벌 비영리 단체가 필리핀의 외딴 지역에서 활동 범위를 넓히려면 어디로 가야할까요?
멸종 위기 조류인 보라색 담비는 북미와 남미를 어떻게 이동하나요? 그 새들은 보호 지역으로 이동하나요?
일본의 어느 지역을 추가적으로 내진 설계를 해야할까요?
캘리포니아의 어느 스타벅스 매장이 다음 스타벅스 리저브 로스터리 매장으로 유력한 후보지인가요?
뉴욕에는 자동차 충돌 사고에 대응할 수 있는 충분한 병원이 있나요? 뉴욕에서 의료 서비스 제공에 공백이 있는 지역은 어디일까요?
또한, 보스턴의 범죄를 시각화하고, 가나의 의료 시설을 조사하고, 유럽의 최고 대학을 탐색하고, 미국의 독성 화학물질 방출을 추적할 수 있습니다.
이 첫번째 튜토리얼에서는 이 강좌를 완료하는데 필요한 전제 조건을 빠르게 다룹니다. 더 깊이 잇는 복습을 원하신다면, Pandas 강좌를 추천합니다.
이제, 첫번째 지리공간 데이터 집합을 시각화해보겠습니다.
데이터 불러오기
첫번째 단계는 지리공간 데이터를 불러오는 것입니다. 이를 위해 GeoPandas 라이브러리를 사용하겠습니다.
지리공간 파일 형식에는 shapefile, GeoJSON, KML, GPKG 등의 다양한 형식이 있습니다. 이 강좌에서는 그 차이점에 대해서는 다루지 않겠지만, 중요한 것들은 다음과 같습니다. - shapefile이 가장 흔한 파일 형식입니다. - 지리공간 파일의 모든 형식은 gpd.read_file() 함수를 이용하여 빠르게 불러올 수 있습니다.
다음 셀은 뉴욕 주 환경보전국에서 관리하는 숲, 야생지대 및 기타 토지에 대한 정보가 포함된 shapefile을 불러옵니다.
# Read in the datafull_data = gpd.read_file("./geospatial-learn-course-data/DEC_lands/DEC_lands/DEC_lands.shp")# View the first five rows of the datafull_data.head()
OBJECTID
CATEGORY
UNIT
FACILITY
CLASS
UMP
DESCRIPTIO
REGION
COUNTY
URL
SOURCE
UPDATE_
OFFICE
ACRES
LANDS_UID
GREENCERT
SHAPE_AREA
SHAPE_LEN
geometry
0
1
FOR PRES DET PAR
CFP
HANCOCK FP DETACHED PARCEL
WILD FOREST
None
DELAWARE COUNTY DETACHED PARCEL
4
DELAWARE
http://www.dec.ny.gov/
DELAWARE RPP
5/12
STAMFORD
738.620192
103
N
2.990365e+06
7927.662385
POLYGON ((486093.245 4635308.586, 486787.235 4...
1
2
FOR PRES DET PAR
CFP
HANCOCK FP DETACHED PARCEL
WILD FOREST
None
DELAWARE COUNTY DETACHED PARCEL
4
DELAWARE
http://www.dec.ny.gov/
DELAWARE RPP
5/12
STAMFORD
282.553140
1218
N
1.143940e+06
4776.375600
POLYGON ((491931.514 4637416.256, 491305.424 4...
2
3
FOR PRES DET PAR
CFP
HANCOCK FP DETACHED PARCEL
WILD FOREST
None
DELAWARE COUNTY DETACHED PARCEL
4
DELAWARE
http://www.dec.ny.gov/
DELAWARE RPP
5/12
STAMFORD
234.291262
1780
N
9.485476e+05
5783.070364
POLYGON ((486000.287 4635834.453, 485007.550 4...
3
4
FOR PRES DET PAR
CFP
GREENE COUNTY FP DETACHED PARCEL
WILD FOREST
None
None
4
GREENE
http://www.dec.ny.gov/
GREENE RPP
5/12
STAMFORD
450.106464
2060
N
1.822293e+06
7021.644833
POLYGON ((541716.775 4675243.268, 541217.579 4...
4
6
FOREST PRESERVE
AFP
SARANAC LAKES WILD FOREST
WILD FOREST
SARANAC LAKES
None
5
ESSEX
http://www.dec.ny.gov/lands/22593.html
DECRP, ESSEX RPP
12/96
RAY BROOK
69.702387
1517
N
2.821959e+05
2663.909932
POLYGON ((583896.043 4909643.187, 583891.200 4...
CLASS 열에서 볼 수 있듯이, 처음 5개 행은 각각 다른 숲에 해당합니다.
이 튜토리얼의 나머지 부분에서는 이 데이터를 사용하여 주말 캠핑 여행을 계획하는 시나리오를 고려해보겠습니다.
온라인에서 크라우드 소싱된 리뷰에 의존하는 대신 자신만의 지도를 만들어봅시다.
이렇게 하면 특정 관심사에 맞게 여행을 조정할 수 있습니다.
전제 조건
데이터의 처음 5개의 행을 보기 위해 head() 함수를 사용했습니다. 이 함수는 Pandas DataFrame에서도 사용된다는 것을 알고 있을 것입니다.
사실, 모든 명령어는 Pandas DataFrame에서도 함께 작동합니다.
이는, 데이터가 (Pandas) DataFrame의 모든 기능을 갖춘 (GeoPandas) GeoDataFrame 객체로 로드되었기 때문입니다.
type(full_data)
geopandas.geodataframe.GeoDataFrame
예를 들어, 모든 열을 사용하지 않으려는 경우 열의 하위 집합을 선택할 수 있습니다. (데이터를 선택하는 다른 방법을 복습하려면, Pandas 강좌의 튜토리얼을 확인하세요.)
data = full_data.loc[:, ["CLASS", "COUNTY", "geometry"]].copy()
value_counts() 함수를 사용하여 다양한 토지 유형 목록과 함께 데이터 집합에 나타나는 횟수를 확인합니다. (이 함수 (관련 함수)를 복습하려면, Pandas 강좌의 튜토리얼을 확인하세요.)
# 각 유형별로 토지의 개수는 몇개인가요?data['CLASS'].value_counts()
이 열에는 다양한 데이터 유형이 포함될 수 있지만, 일반적으로 Point, LineString, or Polygon 입니다.
이 데이터의 geometry 열에는 2,983개의 서로 다른 Polygon 객체가 포함되어 있고, 각 객체는 plot에서 서로 다른 모양에 해당합니다.
아래 셀을 통해 우리는 캠프장 위치 (Point), 도보 경로 (LineString), 지역 경계 (Polygon)를 포함한 3개의 GeoDataFrame을 더 생성합니다.
# 뉴욕 주의 캠프장 위치 (Point)POI_data = gpd.read_file('./geospatial-learn-course-data/DEC_pointsinterest/DEC_pointsinterest/Decptsofinterest.shp')campsites = POI_data.loc[POI_data['ASSET'] =='PRIMITIVE CAMPSITE'].copy()# 뉴욕 주의 도보 경로 (LineString)roads_trails = gpd.read_file('./geospatial-learn-course-data/DEC_roadstrails/DEC_roadstrails/Decroadstrails.shp')trails = roads_trails.loc[roads_trails['ASSET'] =='FOOT TRAIL'].copy()# 뉴욕 주의 지역 경계 (Polygon)counties = gpd.read_file('./geospatial-learn-course-data/NY_county_boundaries/NY_county_boundaries/NY_county_boundaries.shp')
다음으로, 4개의 GeoDataFrame을 모두 사용하여 지도를 만들어봅시다.
plot() 함수는 지도를 커스터마이징하는데 필요로 하는 데 사용할 수 있는 몇가지 파라미터(선택사항)를 입력으로 받습니다.
가장 중요한 것은 ax에 값을 설정하면, 모든 정보가 동일한 지도에 그려진다는 것입니다.
# Base Map을 지역 경계로 설정ax = counties.plot(figsize = (10, 10), color ='none', edgecolor ='gainsboro', zorder =3)# 토지와 캠프장 위치, 도보 경로를 Base Map에 추가wild_lands.plot(color ='lightgreen', ax = ax)campsites.plot(color ='maroon', markersize =2, ax = ax)trails.plot(color ='black', markersize =1, ax = ax)
<Axes: >
북동부 지역 이 캠핑 여행지에 좋은 선택지가 될 것 같습니다!
실습
처음엔 복잡하게 느껴지겠지만, 이미 중요한 분석을 수행할 수 있을만큼 충분히 배웠을 것입니다.
이 강좌에서 만드는 지도는 지구 표면을 2차원으로 묘사합니다. 하지만 아시다시피 지구는 실제로 3차원입니다. 따라서, 지구를 평면으로 렌더링하려면 Map Projection(지도 투영) 이라는 방법을 사용해야합니다.
Map Projection은 100% 정확할 수 없습니다. 각 투영법은 지구 표면을 중요한 속성은 유지하지만, 어떤 식으로든 왜곡합니다.
예를 들어,
equal-area projections(등면적 투영, ‘Lambert Cylindrical Equal Area’, or ‘Africa Albers Equal Area Conic’ 등)은 면적을 보존합니다.
국가나 도시의 면적을 계산하려는 경우 이 투영법을 선택하는 것이 좋습니다.
equidstant projections(등거리 투영, ‘Azimuthal Equidistant projection’)은 거리를 보존합니다.
비행 거리를 계산할 때 좋은 투영법입니다.
투영된 점이 지구상의 실제 위치와 어떻게 일치하는지 보여주기 위해 좌표 참조 체계 (CRS) 를 사용합니다. 이 튜토리얼에서는 좌표 참조계에 대해 자세히 알아보고 GeoPandas에서 좌표 참조계를 사용하는 방법을 알아봅시다.
CRS 설정
shapefile에서 GeoDataFrame을 만들 때, CRS는 이미 로드되어 있습니다.
# 가나의 지역이 포함된 GeoDataFrame 불러오기regions = gpd.read_file('./geospatial-learn-course-data/ghana/ghana/Regions/Map_of_Regions_in_Ghana.shp')print(regions.crs)
epsg:32630
이를 어떻게 해석해야하나요?
좌표 참조 시스템은 유럽 석유 측량 그룹 (EPSG) 코드에 의해 참조됩니다.
이 GeoDataFrame은 일반적으로 ‘Mercator (메르카토르)’ 투영법이라고 더 많이 불리는 EPSG 32630을 사용합니다. 이 투영법은 각도를 보존하고 (해상 항해에 유용) 면적을 약간 왜곡합니다.
그러나, CSV 파일에서 GeoDataFrame을 만들 때, CRS를 설정해야합니다. EPSG 4326은 위도 및 경도 좌표에 해당합니다.
# 가나의 의료 시설로 DataFrame을 만들기facilities_df = pd.read_csv('./geospatial-learn-course-data/ghana/ghana/health_facilities.csv')# DataFrame을 GeoDataFrame으로 변환하기facilities = gpd.GeoDataFrame(facilities_df, geometry = gpd.points_from_xy(facilities_df['Longitude'], facilities_df['Latitude']))# 좌표 참조 시스템 (CRS) EPSG 4326으로 설정facilities.crs = {'init': 'epsg:4326'}# GeoDataFrame의 첫 5개 행 출력facilities.head()
Region
District
FacilityName
Type
Town
Ownership
Latitude
Longitude
geometry
0
Ashanti
Offinso North
A.M.E Zion Clinic
Clinic
Afrancho
CHAG
7.40801
-1.96317
POINT (-1.96317 7.40801)
1
Ashanti
Bekwai Municipal
Abenkyiman Clinic
Clinic
Anwiankwanta
Private
6.46312
-1.58592
POINT (-1.58592 6.46312)
2
Ashanti
Adansi North
Aboabo Health Centre
Health Centre
Aboabo No 2
Government
6.22393
-1.34982
POINT (-1.34982 6.22393)
3
Ashanti
Afigya-Kwabre
Aboabogya Health Centre
Health Centre
Aboabogya
Government
6.84177
-1.61098
POINT (-1.61098 6.84177)
4
Ashanti
Kwabre
Aboaso Health Centre
Health Centre
Aboaso
Government
6.84177
-1.61098
POINT (-1.61098 6.84177)
위의 셀 CSV 파일에서 GeoDataFrame을 만들려면 Pandas와 GeoPandas를 모두 사용해야했습니다.
먼저 위도 및 경도 좌표가 포함된 열을 포함하는 DataFrame 만듭니다.
이를 GeoDataFrame으로 변환하기 위해 gpd.GeoDataFrame()을 사용합니다.
gpd.points_from_xy() 함수는 위도 및 경도 열에서 Point 개체를 생성합니다.
재투영 (Re-projecting)
재투영은 CRS를 변경하는 과정을 말합니다. 이 작업은 GeoPandas에서 to_crs() 함수를 사용하여 수행됩니다.
여러 개의 GeoDataFrame을 plot할 때 모두 동일한 CRS를 사용하는 것이 중요합니다. 아래 셀에서는 facilities GeoDataFrame의 CRS를 regions의 CRS와 일치하도록 변경한 후 plot합니다.
# 지도 만들기ax = regions.plot(figsize = (8, 8), color ='whitesmoke', linestyle =':', edgecolor ='black')facilities.to_crs(epsg =32630).plot(markersize =1, ax = ax)
<Axes: >
to_crs() 함수는 geometry 열만 수정하고 다른 모든 열은 그대로 유지합니다.
# `Latitude`와 `Longitude` 열은 변경 Xfacilities.to_crs(epsg =32630).head()
Region
District
FacilityName
Type
Town
Ownership
Latitude
Longitude
geometry
0
Ashanti
Offinso North
A.M.E Zion Clinic
Clinic
Afrancho
CHAG
7.40801
-1.96317
POINT (614422.662 818986.851)
1
Ashanti
Bekwai Municipal
Abenkyiman Clinic
Clinic
Anwiankwanta
Private
6.46312
-1.58592
POINT (656373.863 714616.547)
2
Ashanti
Adansi North
Aboabo Health Centre
Health Centre
Aboabo No 2
Government
6.22393
-1.34982
POINT (682573.395 688243.477)
3
Ashanti
Afigya-Kwabre
Aboabogya Health Centre
Health Centre
Aboabogya
Government
6.84177
-1.61098
POINT (653484.490 756478.812)
4
Ashanti
Kwabre
Aboaso Health Centre
Health Centre
Aboaso
Government
6.84177
-1.61098
POINT (653484.490 756478.812)
GeoPandas에서 EPSG 코드를 사용할 수 없는 경우, CRS의 ’proj4 string’을 사용하여 CRS를 변경할 수 있습니다.
folium.Marker() 함수로 지도에 Marker를 추가합니다. 아래의 각 Marker는 서로 다른 강도에 해당합니다.
# 지도 생성m_2 = folium.Map( location = [42.32, -71.0589], tiles ='cartodbpositron', zoom_start =13)# 지도 + Markerfor idx, row in daytime_robberies.iterrows() : Marker([row['Lat'], row['Long']]).add_to(m_2)# 지도 표현m_2
Make this Notebook Trusted to load map: File -> Trust Notebook
folium.plugins.MarkerCluster()
folium.plugins.MarkerCluster()를 사용하면 맵을 깔끔하게 정리할 수 있습니다. 각 Marker는 MarkerCluster 객체에 추가됩니다.
# 지도 생성m_3 = folium.Map( location = [42.32, -71.0589], tiles ='cartodbpositron', zoom_start =13)# 지도 + MarkerClustermc = MarkerCluster()for idx, row in daytime_robberies.iterrows() :ifnot math.isnan(row['Long']) andnot math.isnan(row['Lat']): mc.add_child(Marker([row['Lat'], row['Long']]))m_3.add_child(mc)# Display the mapm_3
Make this Notebook Trusted to load map: File -> Trust Notebook
버블맵 (Bubble Maps)
버블 맵는 Marker 대신 원을 사용합니다. 각 원의 크기와 색상을 변경하여 위치와 다른 두 변수 사이의 관계를 표시할 수도 있습니다.
folium.Circle()을 사용하여 원을 반복적으로 추가하여 버블 지도를 만듭니다.
아래 셀에서 9~12시에 발생한 강도는 녹색으로 표시되고, 13~17시에 발생한 강도는 빨간색으로 표시됩니다.
# 지도 생성m_4 = folium.Map(location = [42.32, -71.0589], tiles ='cartodbpositron', zoom_start =13)def color_producer(val) :if val <=12 :return'forestgreen'else :return'darkred'# 지도 + 버블맵for i inrange(0, len(daytime_robberies)) : Circle( location = [daytime_robberies.iloc[i]['Lat'], daytime_robberies.iloc[i]['Long']], radius =20, color = color_producer(daytime_robberies.iloc[i]['HOUR'])).add_to(m_4)# 지도 표현m_4
Make this Notebook Trusted to load map: File -> Trust Notebook
folium.Circle()은 여러 인수를 받습니다.
location은 원의 중심을 위도와 경도로 포함하는 목록입니다.
radius는 원의 반지름을 설정합니다.
기존 버블맵에서는 각 원의 반지름이 달라질 수 있습니다. 각 원의 색상을 변경하는데 사용되는 color_producer() 함수와 유사한 함수를 정의하여 이를 구현할 수 있습니다.
color는 각 원의 색상을 설정합니다.
color_producer() 함수는 강도의 위치에 대한 시간의 효과를 시각화하는데 사용됩니다.
히트맵 (Heatmaps)
히트맵을 만드려면 folium.plugins.HeatMap()을 사용합니다. 이는 도시 내 여러 지역의 범죄 밀도를 보여주며, 빨간색 영역은 상대적으로 범죄 발생이 더 많습니다.
대도시에서 예상할 수 있듯이, 대부분의 범죄는 도심 근처에서 발생합니다.
# 지도 생성m_5 = folium.Map(location = [42.32, -71.0589], tiles ='cartodbpositron', zoom_start =12)# 지도 + 히트맵HeatMap(data = crimes[['Lat', 'Long']], radius =10).add_to(m_5)# 지도 표현m_5
Make this Notebook Trusted to load map: File -> Trust Notebook
folium.plugins.HeatMap()은 여러 인수를 받습니다.
data는 시각화하고자 하는 위치를 포함하는 DataFrame입니다.
radius는 히트맵의 부드러움 정도를 조정합니다. 값이 클수록 히트맵의 부드러움 정도가 커집니다.(즉, 간격이 줄어듭니다.)
단계구분도 (Choropleth Maps)
경찰 관할 구역별로 범죄가 어떻게 다른지 이해하기 위해 단계구분도를 만들어 보겠습니다.
첫번째 단계로, 각 구역에서 서로 다른 행이 할당되고 geometry 열에 지리적 경계가 포함된 GeoDataFrame을 만듭니다.
# 보스턴 경찰 관할 구역 지리적 경계가 포함된 GeoDataFramedistricts_full = gpd.read_file('./geospatial-learn-course-data/Police_Districts/Police_Districts/Police_Districts.shp')districts = districts_full[['DISTRICT', 'geometry']].set_index('DISTRICT')districts.head()
geometry
DISTRICT
A15
MULTIPOLYGON (((-71.07416 42.39051, -71.07415 ...
A7
MULTIPOLYGON (((-70.99644 42.39557, -70.99644 ...
A1
POLYGON ((-71.05200 42.36884, -71.05169 42.368...
C6
POLYGON ((-71.04406 42.35403, -71.04412 42.353...
D4
POLYGON ((-71.07416 42.35724, -71.07359 42.357...
각 구역별 범죄 발생 건수를 보여주는 plot_dict라는 Pandas Series를 만듭니다.
# 각 구역별 범죄 발생 건수plot_dict = crimes['DISTRICT'].value_counts()plot_dict.head()
data는 각 지리적 영역의 색상을 지정하는데 사용할 값이 포함된 Pandas Series입니다.
key_on은 항상 feature.id로 설정됩니다.
이는 geo_data에 사용되는 GeoDataFrame과 data에 사용되는 Pandas Series가 동일한 인덱스를 가지고 있다는 사실을 나타냅니다. 자세한 내용을 이해하려면 GeoJSON Feature Collection의 구조를 좀 더 자세히 살펴봐야합니다. (‘feature’ 키에 해당하는 값은 list이고, 각 항목은 ‘id’ 키가 포함된 dictonary입니다.)
각 국가의 예상 인구와 국내총생산(GDP)을 포함하는 DataFrame europe_stats와 결합하겠습니다.
europe_stats.head()
name
pop_est
gdp_md_est
0
Russia
144373535.0
1699876
1
Norway
5347896.0
403336
2
France
67059887.0
2715518
3
Sweden
10285453.0
530883
4
Belarus
9466856.0
63080
아래 셀에서 특성 조인을 수행합니다. on 인수는 europe_boundaries의 행을 europe_stats의 행과 일치시키는데 사용되는 열 이름으로 설정됩니다.
# 특성 조인을 사용하여 유럽 국가에 대한 데이터 병합europe = europe_boundaries.merge(europe_stats, on ='name')europe.head()
name
geometry
pop_est
gdp_md_est
0
Russia
MULTIPOLYGON (((180.00000 71.51571, 180.00000 ...
144373535.0
1699876
1
Norway
MULTIPOLYGON (((15.14282 79.67431, 15.52255 80...
5347896.0
403336
2
France
MULTIPOLYGON (((-51.65780 4.15623, -52.24934 3...
67059887.0
2715518
3
Sweden
POLYGON ((11.02737 58.85615, 11.46827 59.43239...
10285453.0
530883
4
Belarus
POLYGON ((28.17671 56.16913, 29.22951 55.91834...
9466856.0
63080
공간 조인 (Spatial Join)
또 다른 조인 유형은 공간 조인(Spatial Join)입니다. 공간 조인을 사용하면 geometry 열에 있는 개체 간의 공간 관계를 기반으로 GeoDataFrame을 결합합니다. 예를 들어, 유럽 대학의 지오코딩된 주소가 포함된 GeoDataFrame universities가 이미 있습니다.
그런 다음 공간 조인을 사용하여 각 대학을 해당 국가에 일치시킬 수 있습니다. 이 작업은 gpd.sjoin()을 사용하여 수행합니다.
# 공간 조인을 사용하여 대학을 유럽의 국가와 일치european_universities = gpd.sjoin(universities, europe)# 결과 검토print('We located {} universities.'.format(len(universities)))print('Only {} of the universities were located in Europe (in {} different countries).'.format(len(european_universities), len(european_universities.name.unique())))european_universities.head()
We located 91 universities.
Only 86 of the universities were located in Europe (in 14 different countries).
Name
Latitude
Longitude
geometry
index_right
name
pop_est
gdp_md_est
0
University of Oxford
51.759037
-1.252430
POINT (-1.25243 51.75904)
28
United Kingdom
66834405.0
2829108
1
University of Cambridge
52.200623
0.110474
POINT (0.11047 52.20062)
28
United Kingdom
66834405.0
2829108
2
Imperial College London
51.498959
-0.175641
POINT (-0.17564 51.49896)
28
United Kingdom
66834405.0
2829108
4
UCL
51.521785
-0.135151
POINT (-0.13515 51.52179)
28
United Kingdom
66834405.0
2829108
5
London School of Economics and Political Science
51.514211
-0.116808
POINT (-0.11681 51.51421)
28
United Kingdom
66834405.0
2829108
위의 공간 조인은 두 GeoDataFrame의 geometry 열을 살펴봅니다. universities GeoDataFrame의 Point 객체가 europe DataFrame의 Polygon 객체와 교차하는 경우, 해당 행이 결합되어 european_universities DataFrame의 단일 행으로 추가됩니다. 그렇지 않으면, 일치하는 대학이 없는 국가(및 일치하는 국가가 없는 대학)는 결과에서 제거됩니다.
gpd.sjoin() 함수는 how 및 op 인수를 통해 다양한 조인 유형에 맞게 지정할 수 있습니다. 예를 들어, how = left(또는 how = right)를 설정하여 SQL Left (또는 Right) 조인과 동등한 작업을 수행할 수 있습니다. 이 강좌에서는 자세히 설명하지 않겠지만, 이 문서에서 자세한 내용을 확인할 수 있습니다.
서로 다른 두 GeoDataFrame의 Point 간의 거리를 측정하려면, 먼저 두 GeoDataFrame이 동일한 좌표 참조 시스템(CRS)을 사용하는지 확인해야 합니다.
다행히도 여기서는 둘 다 EPSG:2272를 사용합니다.
print(stations.crs)print(releases.crs)
epsg:2272
epsg:2272
또한, CRS가 어떤 단위(meter, feet 또는 다른 단위)를 사용하는지 확인합니다. 이 경우 EPSG:2272는 feet 단위를 사용합니다. (원한다면 여기에서 확인할 수 있습니다.)
GeoPandas에서 거리를 계산하는 것은 비교적 간단합니다. 아래 셀은 recent_release의 비교적 최근 릴리즈 사건과 stations GeoDataFrame의 모든 스테이션 사이의 거리(피트)를 계산합니다.
# 특정 릴리즈 사건 1개 선택recent_release = releases.iloc[360]# 릴리즈에서 각 스테이션까지의 거리 측정distances = stations['geometry'].distance(recent_release['geometry'])distances
각 Polygon을 지도에 그리기 위해 folium.GeoJson()을 사용합니다. folium에는 위도와 경도의 좌표가 필요하므로 시각화 전에 CRS를 EPSG:4326으로 변환해야합니다.
# 최근 릴리즈 사건과 모니터링 스테이션을 포함하여 지도를 생성m = folium.Map(location = [39.9526, -75.1652], zoom_start =11)HeatMap(data = releases[['LATITUDE', 'LONGITUDE']], radius =15).add_to(m)for idx, row in stations.iterrows() : Marker([row['LATITUDE'], row['LONGITUDE']]).add_to(m)# 지도 + 버퍼GeoJson(two_mile_buffer.to_crs(epsg =4326)).add_to(m)# 지도 표현m
Make this Notebook Trusted to load map: File -> Trust Notebook
이제 모든 모니터링 스테이션에서 2마일 이내에 독성 물질 방출이 발생했는지 테스트하려면 각 Polygon에 대해 12가지 테스트를 실행하여 해당 지점이 포함되어 있는지 개별적으로 확인할 수 있습니다.
하지만 더 효율적인 방법은 먼저 모든 Polygon을 MultiPolygon 객체로 축소하는 것입니다. 이 작업은 unary_union 속성을 사용하여 수행합니다.
# Polygon 그룹을 단일 MultiPolygon 객체로 축소my_union = two_mile_buffer.geometry.unary_unionprint('Type : ', type(my_union))# MultiPolygon 객체 시각화my_union
Type : <class 'shapely.geometry.multipolygon.MultiPolygon'>
MultiPolygon에 Point가 포함되어 있는지 확인하기 위해 contains() 함수를 사용합니다. 튜토리얼 앞부분의 릴리즈 사건을 사용하겠습니다.
가장 가까운 모니터링 스테이션까지의 거리가 약 3,781피트라는 것을 알고 있습니다.
하지만 모든 릴리즈가 대기 모니터링 스테이션에서 2마일 이내에 발생한 것은 아닙니다!
# 가장 가까운 스테이션이 2마일 미만인 경우print(my_union.contains(releases.iloc[360].geometry))# 가장 가까운 스테이션이 2마일 이상 떨어진 경우print(my_union.contains(releases.iloc[358].geometry))